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New evaluation of the relativistic effects in the double production of S*-wave char- 
monium states is performed on the basis of perturbative QCD and the relativistic 
quark model. The main improvement consists in the exact account of properties 
of the relativistic meson wave functions. For the gluon and quark propagators en- 
£NJ ■ tering the production vertex function we use a truncated expansion in the ratio of 

(T) ■ the relative quark momenta to the center-of-mass energy ^/s up to the second or- 

der. The exact relativistic treatment of the wave functions makes all such second 

00 

\ order terms convergent, thus allowing the reliable calculation of their contributions 

to the production cross section. Compared to the nonrelativistic calculation we ob- 
tain a significant increase of the cross sections for the S'-wave double charmonium 
production. This brings new theoretical results in good agreement with the available 
experimental data. 
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The production processes of mesons and baryons containing heavy b and c quarks in 
different reactions are under intensive study at present [H, 0, U 0] • The experimental in- 
vestigation of the double charmonium production in e + e~ annihilation by BaBar and Belle 
Collaborations revealed a discrepancy between the measured cross sections and theoreti- 
cal results obtained in the nonrelativistic approximation in QCD @, 0] . Various efforts 
have been undertaken to improve the theoretical calculations. They include the evaluation 
of radiative corrections of order a s and the investigation of relativistic effects connected 
with the relative motion of the heavy quarks forming the vector and pseudoscalar quarkonia 
@, S, 0, 0, El El, H H E EH, EJ H HH. Asa result, the difference between theory 



and experiment for the value of the center-of-mass energy a/s = 10.6 GeV was essentially 
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FIG. 1: The production amplitude of a pair of charmonium states (V denotes the vector meson 
and V the pseudoscalar meson) in e + e~ annihilation. The wave line shows the virtual photon and 
the dashed line corresponds to the gluon. V is the production vertex function. 



decreased [3, E, 13, 15]. Moreover, the new theoretical analysis carried out in Refs. 

Urn 

shows that the inclusion of order a s and relativistic corrections decreases the discrepancy 
between theory and experiment at the present level of precision. But despite this fact there 
exists the frequently debated question connected with the calculation of the relativistic cor- 
rections in the production cross section. It is related to the determination of the specific 
parameter (p 2 )=J p 2 ^' V (p)dp/ (2n) 3 emerging after the expansion of all quantities in the 
production amplitude in the relative quark momenta p and q 15, 20, 21, 22|, where are 
the vector and pseudoscalar charmonium wave functions in the rest frame. The divergence 
of this integral required the use of a regularization procedure (dimensional regular izat ion is 
commonly used) which led to a definite uncertainty of the evaluation. Moreover, the large 
value of the relativistic contribution obtained in the previous studies 0, lEf evidently rises 
a question about the convergence of the expansion in the heavy quark velocity. In this 
letter we propose an alternative approach to the calculation of relativistic effects based on 



the relativistic quark model [23J, |24], |25|, |26|, |27| and perturbative QCD. It uses a truncated 



expansion in relative momenta p and q and thus avoids divergent integrals in the relativistic 
contribution of the second order. 

Define the four momenta of the produced c, c quarks forming the vector and pseudoscalar 
charmonia in terms of total momenta P(Q) and relative momenta p(q) as follows: 



Pl,2 



1 



-P±p, (p-P) = 0; g 1)2 = -Q±g, (q ■ Q) = 



(1) 



where p = Lp(0,p), q = Lp(Q, q) are the four-momenta obtained from the rest frame 
four-momenta (0, p) and (0, q) by the Lorentz transformation to the system moving with 
the momenta P, Q. Then the production amplitude of the S- wave vector and pseudoscalar 



charmonium states, shown in Fig.l, can be presented in the form [15|, |28|, |29 
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where a s is the QCD coupling constant, a is the fine structure constant, Q c is the c quark 
electric charge. The relativistic wave functions of the bound quarks ^ v > v accounting for the 
transformation from the rest frame to the moving one with four momenta P, Q are 
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where t>i = P/My, t> 2 = Q/M-p] i is the polarization vector of the vector charmonium; 
e(p) = Vp 2 + m2 and m is the c quark mass. The vertex function T v {jp, P; q, Q) at leading 
order in a, can be written as a sum of four contributions: 



r»(p,P;q,Q) = 1 , 



(f — q~i + m) 



- 1/3 D^{k 2 )+ lp - 



r — q\Y — m A + le 
(r — q 2 ) z — m z + xe i 



(pi — r + m) 
— pi) 2 — m 2 + ie 
(r - p 2 + m) 



r — p 2 ) — m 2 + ie 
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where the gluon momenta are k\ = p\ + qi, k 2 = p 2 + q 2 and r 2 = s = (P + Q) 2 = (p_ +p + ) 2 , 
P-, p + are four momenta of the electron and positron. The dependence on the relative 
momenta of c-quarks is present both in the gluon propagator D^ v {k) and quark propagators 
as well as in the relativistic wave functions. One of the main technical difficulties in calcu- 
lating the production amplitude ([2]) consists in performing angular integrations, since both 
gluon and quark propagators in the vertex function (jSj) contain angles in the denominators. 
Therefore we expand these propagators in the relative momenta. Such expansion leads to 
the vertex function containing angles only in numerators and, thus, the angular integrations 
can be easily performed. 

The inverse denominators of quark propagators expanded in the ratio of the relative quark 
momenta p, q to the energy a/s up to the second order can be expressed as follows: 
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where the factors Z\ and Z 2 differ only due to the bound state corrections: 
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Corresponding expansions of the gluon propagators in Eq.(5) with the account of terms of 
order 0(p 2 /s,q 2 / s) are (Z = s/4): 
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We expanded the gluon and quark propagators in the ratio of the relative quark momenta 
to the center-of-mass energy ^Js up to the second order terms in the production vertex 
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function (5) but preserved all relativistic factors entering the denominators of the relativistic 
wave functions (3), (4). This provides the convergence of the resulting momentum integrals. 
Then keeping the terms of second and fourth order in both variables p and q in the numerator 
of Eq.(2) from the relativistic wave functions (3)-(4) and second order from the expansions 
of the quark and gluon propagators, we perform the angular averaging taking into account 
Eq.(l) and using the following relation: 
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Then we can write the total production amplitude M. in the form: 
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where T13 originates from the sum of the first and third terms in the vertex function (5 



a3 



and T 2 4 from the sum of the second and fourth terms. First, using the Form package 
we presented T 13 and T 24 as a series over the factors Zx, Z 2) u = M v j{M v + My), n = 
m/(M P + M v ), c(p) = [2m/(e(p)+m)-l] = -p 2 /(e(p) + m) 2 , c(q) = [2m/(e(q) +m) - 1] = 
— q 2 / (e(g) +m) 2 . The resulting expressions are cumbersome so we omit them here. 1 Then we 
performed their simplification by neglecting the bound state corrections in the denominators 
Z\ and Z2 (8). This can be done because the value of \fs at which the experimental data 
were obtained is essentially larger than the quark bound state energy. In this approximation 
which does not influence the accuracy of the calculation (the corresponding error in the cross 
section at the energy yfs = 10^-11 GeV amounts 0.5 %) we have Z\ Z 2 ~ s/2. After such 
approximation the total cross section for the exclusive production of pseudoscalar and vector 
charmonium states in e + e~ annihilation is given by the following analytical expression: 
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where the function T(p, q) can be written as follows: 
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The nonzero values of the coefficients u>m, Pkh &kh 7i,2 are given explicitly in the Appendix A. 

The momentum integrals entering Eq.(12) are convergent and we calculate them nu- 
merically, using the wave functions obtained by the numerical solution of the relativistic 



They are available from authors: apm@physik.hu-berlin.de 
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FIG. 2: The cross section in fb of e + e~ annihilation into a pair of S-wave charmonium states with 
the opposite charge parity as a function of the center-of-mass energy ^/s (solid line). The dashed 
line shows the nonrelativistic result without bound state and relativistic corrections. 



quasipotential wave equation [2J, |25|, |26|, |31(. The exact form of the wave functions \l/ (p) 
and ^(q) is extremely important for getting the reliable numerical results. It is sufficient to 
note that the charmonium production cross section a(s) in the nonrelativistic approximation 
contains the factor l^^flWPl^jv/ilC')! 2 - So, small changes of the numerical values ^^(0) 
and ty^ R (0) considerably influence the final result. In the approach based on nonrelativistic 
QCD this problem is closely related to the determination of the color-singlet matrix elements 
for the charmonium. Therefore for our calculations we use the charmonium wave functions 
^) V ' V obtained with the complete nonperturbative treatment of relativistic effects. For this 
purpose we consider the quark-antiquark interaction operator constructed in the relativistic 
quark model in Refs. ^2a 



25. 26 



Thus, in the present study of the production amplitude 
(2) we keep the relativistic corrections of two types. The first type is determined by several 
functions depending on the relative quark momenta p and q arising from the gluon propa- 
gator, the quark propagator and the relativistic meson wave functions. The second type of 
corrections originate from the nonperturbative treatment of the hyperfine interaction in the 
quark-antiquark potential which leads to the different wave functions \1/q (p) and ^(q) for 
the vector and pseudoscalar charmonium states, respectively. In addition, we systematically 
accounted the bound state corrections working with the observed masses of the vector and 
pseudoscalar mesons. The calculated masses of vector and pseudoscalar charmonium states 
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TABLE I: Comparison of the obtained results with previous theoretical predictions and experi- 
mental data. 



State 
H1H2 


a BaBar^ 
BrH 2 ^charged>2 

m m 


^ Belief 
BrH 2 ^charged>2 

(fb) [5] 


[12] 


& NRQCD 

(fb) m 




[7J 


[15] 


Our result 
(fb) 


*(15)t/ c (15) 


17.6 ± 2.8±J;f 


25.6 ±2.8 ±3.4 


26.7 


3.78 


5.5 


7.4 


7.8 


22.2 ±3.6 


*(2S) Vc (lS) 




16.3 ±4.6 ±3.9 


16.3 


1.57 


3.7 


6.1 


6.7 


15.3 ±2.4 


*(1S)7/ C (2S) 


16.4 ± 3.7t|o 


16.5 ±3. ±2.4 


26.6 


1.57 


3.7 


7.6 


7.0 


16.4 ±2.6 


*(25)t/ c (25) 




16.0 ±5.1 ±3.8 


14.5 


0.65 


2.5 


5.3 


5.4 


9.6 ± 1.5 



agree well with experimental values [25|, [32] . Note that all parameters of the model are kept 
fixed from the previous calculations of the meson mass spectra and decay widths 24j, |25|, [27 . 
The masses of the S'-wave charmonium states are: mj/q, = 3.097 GeV, m Vc = 2.980 GeV, 
m^i = 3.686 GeV, m v ' c = 3.637 GeV. The strong coup ling; constant entering the production 
amplitude (2) is taken to be a s =0.21 (see also |7I. Il2j]). 

Numerical results and their comparison with several previous calculations and experimen- 
tal data are presented in Table I. In Refs.[l~7. lH 20] the cross section a[e + + e~ — > J/^f + T] c ] 
was calculated with the values 20.04 fb, 17.5 ± 5.7 fb and 17.6^^7 fb, respectively. The cal- 
culated production cross sections of a pair of S'-wave charmonium states are shown in Fig. 2. 
Our new evaluation of the cross sections in the reaction e + + e~ — > V cg + V cS evidently shows 
that the systematic account of all relativistic effects connected with the bound state wave 
functions, the gluon and quark propagators removes the discrepancy between theory and 
experiment. Numerically, the increase of the cross section a (12) is determined approxi- 
mately by the factor of 2 coming from the relativistic corrections entering in the production 
amplitude (2) (in this part our results agree with the previous calculations in Ref . [l5| ) and 
by another factor of 2 from the relativistic bound state wave functions. In our analysis we 
use the exact expressions (3)-(4) for the relativistic wave functions. Thus we correctly take 
into account all relativistic contributions of orders 0(v 2 ) and 0(v 4 ) since they are deter- 
mined by the convergent momentum integrals due to the presence of the relativistic factors 
in the denominators of expressions (3)-(4). Therefore the resulting theoretical uncertainty is 
connected with the omitted terms of the employed truncated expansions (6), (7), (9) which 
are of order v 2 p 2 /s. Taking into account that the average value of the heavy quark velocity 
squared in the charmonium is < v 2 >= 0.3, we expect that they should not exceed 5-10% in 
the interval of energies y/s = 7-7-11 GeV. We should remind also that our relativistic quark 
model has the phenomenological structure and differs significantly from the approach of non- 
relativistic QCD (NRQCD). Despite the fact that it is based on the quantum field-theoretic 
approach, it contains a number of the phenomenological parameters which we fixed solving 
many tasks in the quarkonium physics. Unfortunately, we can not control the theoretical 
accuracy in the same manner as in NRQCD. We obtained the theoretical predictions for the 
masses and decay rates of different charmonium states with more than one per cent accuracy. 
So, we suppose in this study that there are no additional essential theoretical uncertainties 
in the bound state wave functions connected with the formulation of our model in regions 
of nonrelativistic and relativistic momenta. 

It is important to point out that it is not possible to simply compile the enhancements 
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of the production cross sections originating from our calculation of the relativistic contribu- 
tions and from the one-loop corrections calculated in Ref . [l3| . The latter was done in the 
nonrelativistic limit. Indeed in our model the interaction potential in the relativistic wave 
equation contains the one-loop radiative corrections. Therefore the inclusion of the one- 
loop corrections considered in [131] in our calculation requires their complete recalculation 
using our relativistic wave functions, since we take into account effectively some part of the 
one-loop diagrams connected with the exchange of gluons between heavy quarks in the final 
state 2 . As a result both nonperturbative and partially perturbative contributions are taken 
into account. 

Thus our approach cannot be directly confronted with the one of Ref. [l3| . The radiative 
corrections are the main source of the theoretical uncertainty in our calculations. Indeed, 
available estimates of one-loop corrections in the nonrelativistic limit indicate that they are 
considerable. Taking their values from [l3[ (relativistic factor K — 1.8 to nonrelativistic 
result) we estimate that this part of the theoretical error should not exceed 15%. Therefore 
the total theoretical uncertainty amounts to 16% for the energy region y/s = 10.6 GeV. To 
obtain this estimate we add the above mentioned relativistic and one-loop uncertainties in 
quadrature (as it was done in Ref.[l~9l]). These theoretical errors in the calculated production 
cross section at yfs = 10.6 GeV are shown directly in Table I. There are no additional 
uncertainties related to the choice of m c or any other parameters of the model, since their 
values were fixed from our previous consideration of meson and baryon properties 24|, |25, 



26, 27, 29 



In summary, we presented a systematic treatment of relativistic effects in the double 
charmonium production in e + e~ annihilation. We explicitly separated two different types of 
relativistic contributions to the production amplitudes. The first type includes the relativis- 
tic v/c corrections to the wave functions and their relativistic transformations which were 
for the first time exactly taken into account. The second type includes the relativistic pj \/s 
corrections emerging from the expansion of the quark and gluon propagators. The latter 
corrections were taken into account up to the second order. It is important to note that 
the expansion parameter p/y/s is very small. Contrary to the previous calculations within 
NRQCD all obtained expressions for the relativistic contributions are now expressed through 
converging integrals. Thus no additional uncertainty related to their regularization emerges. 
Therefore we can reliably estimate the uncertainty originating from the neglected higher- 
order relativistic contributions. The calculated values for the production cross sections agree 
well with experimental data. 
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APPENDIX A: THE COEFFICIENTS u^, a ijt Pij , 7i ENTERING IN THE 

PRODUCTION CROSS SECTION 



^01 



Uoo = -18(-1 + 2k - 3u)(u - 1)V, 
6(u - l) 2 [32k 3 + 16k 2 (-5 + u) - 6km 2 + 3(1- 5m) u 2 



ujw = 6m 



(Al) 
(A2) 

(A3) 

(A4) 

(A5) 
(A6) 
(A7) 



27)) 



(A8) 



96k 3 - 34k (u - if + (u- if (5m - 1) - 16/t 2 (1 + 11m) 
wii = -2 [536k 5 + 102k (u - if u 2 + 3 (u - if u 2 (1 + 3m) - 8k 4 (61 + 67m) 
-6k 3 (114 + u (161m - 228)) - 2k 2 (46 + u (-170 + u (175 + 237m))) 
Ul2 = -96km 2 [6k 2 - 2(m - l) 2 - k(1 + 11m)] , 
uj 21 = -96k 2 (m - l) 2 (-5 + 2k + u) , 
7l = 64k 3 [64k 6 (5 + 2m (2m - 5)) + 32k 5 (4 + (m - 1) m (14 + 3m)) + 
+4k 2 (m - l) 2 (146 + m (-584 + u (851 + 6m (21m - 89)))) - 16k 4 x 
(44 + m (m (273 + m (49m - 194)) - 176))-16k 3 (7 + u (u (37 + u (u - 25 + 15m 2 )) 

-(m - l) 2 (194 + m (-1164 + m (2754 + u (-3256 + 3m (676 + 27 (u - 8) «))))) + 
+2k(m - 1) (2 + m (-108 + m (523 + u (-992 + u (874 + 3m (9m - 106)))))) 

72 = -128k 4 (1 - 2m) 2 [-31 + 32k 5 + 2k(1 - 3m) 2 (1 + m) 2 - 16k 4 (1 + 3m) - 
-16k 3 (1 + m (5m - 2)) + 8k 2 (7 + u (-17 + u (11 + 15m))) + 
+m (151 - m (286 + 3m (-98 + u (55 + 9m))))] , 

a 01 = 96k 2 (m - l) 2 [-17 + 2k(1 - 2m) 2 + 16k 4 (-1 + u) - 
-8k 2 (m - 1)(4 + m(5m - 8)) + m (81 + u (-144 + u (l24 - 57m + 9m 2 ) 
a 10 = 48k 2 m 2 [-1 + 32k 5 + 16k 4 (1 - 5m) - 16k 3 (1 + m(5m - 2)) + 
+8k 2 (5m - 1)(1 + m(5m - 2)) + 6k (1 + u (3m - 2) (2 + u (2 + u))) - 
-m (7 + m (-34 + m (2 + 51m + 45m 2 )))] , 

a u = k 2 [600 - 8k (-1695 + 2k (-45 + k (1527 + 2k (11 + k (-659 + 2k (67k - 31)))))) - 

(All) 

-4328m + 16k (-5085 + k (-769 + 2k (3054 + k (-351 + 2k (97k - 659))))) u+ 
+4 (3235 + k (49319 + 4k (2783 + k (-9443 + k (1071 + 1681k))))) m 2 - 
-4(4987 + 4k(15419 + k(3801 + 667k(k - 10))))m 3 + 2(7933 + k(83683+ 

+2 (7829 - 8265k) k))m 4 - 2 (3711 + 29142k + 3770k 2 ) m 5 + (4765 + 9851k)m 6 + 983m 7 ] , 

p 01 = 96k 2 (m - l) 2 [l6 + 12k 2 (m - 1) + m(-40 + (33 - 13m)m) - 2k(2 + (m - 4)u)l , (A12) 



(A9) 
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p 10 = 24ku 2 [48k 4 + 16k 3 (1 - 7m) + k 2 (-4 + 8u - 76u 
-{u- l) 2 (-7 + u{U + 3m)) + 2k (-1 + u (-1 + u + 73m 2 



(A13) 



p u = -8k [528K b - 128k 5 (3 + 5m) + 3(m - 1) V(-7 + u(14 + 3m)) - 12/t 4 (67 + 2m(56m - 67)) + 

(A14) 

+4k 3 (97 + m (-259 + 6m(39 + 9m))) + k 2 (8U + m(-3256 + m(5017 + m(-3522 + 1159m)))) + 
+k (70 + m (-382 + m (969 + u (-1323 + 2m(469 + 80)))))] , 

P21 — —POI, Pl2 — —Rio, °21 — — Coi) °\2 — — °"10- 
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